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LONG-TERM  GOALS 

The  overall  goal  of  this  project  is  to  design  adaptive  array  processing  algorithms  that  have  good 
transient  performance,  are  robust  to  mismatch,  work  with  low  sample  support,  and  incorporate 
waveguide  physics  in  a  meaningful  and  robust  way. 

OBJECTIVES 

Adaptive  array  processing  algorithms  facilitate  the  detection  and  localization  of  quiet  sources  by  nulling 
out  noise  and  interference.  These  algorithms,  which  use  the  incoming  data  to  design  optimal  weight 
vectors,  provide  substantial  gains  in  performance  over  non-adaptive  techniques.  When  the  input  is 
non-stationary,  however,  the  performance  of  adaptive  processors  may  be  significantly  degraded  due  to 
low  sample  support.  Since  ocean  acoustic  signals  are  often  non-stationary  due  to  a  number  of  factors, 
e.g .,  source  motion,  receiver  motion,  and  environmental  fluctuations,  it  is  crucial  to  have  algorithms  that 
work  in  rapidly  changing  environments.  Previous  work  has  typically  focused  on  the  asymptotic 
performance  of  adaptive  processors,  and  much  less  attention  has  been  given  to  the  case  where  the  input 
signals  change  faster  than  the  processor  can  reach  the  asymptotic  limit.  This  project  focuses  on  the 
general  problem  of  array  processing  in  realistic  non-stationary  ocean  environments.  It  has  several 
specific  objectives:  1)  To  characterize  the  performance  of  existing  adaptive  algorithms,  focusing 
specifically  on  the  transient  behavior  and  convergence  rate,  rather  than  asymptotic  performance;  2)  To 
design  new  adaptive  algorithms  that  are  robust  and  have  the  rapid  convergence  that  is  required  in 
non-stationary  environments;  3)  To  explore  methods  of  incorporating  propagation  physics  into  adaptive 
processing  algorithms  without  sacrificing  robustness. 

APPROACH 

To  develop  improved  array  processing  algorithms  for  non-stationary  environments,  this  project  is 
investigating  adaptive  beamforming  algorithms  that  can  operate  on  a  single-snapshot  or  a  very  limited 
number  of  snapshots.  Spectral  estimation  methods,  such  as  the  multitaper  spectral  estimator  developed 
by  Thomson  [1],  require  only  a  single  snapshot  to  produce  an  estimate.  This  project  is  exploring  how  to 
adapt  the  multitaper  approach  to  the  spatial  spectral  estimation  problem. 
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The  signal  processing  algorithms  developed  in  this  project  are  being  tested  on  both  simulated  and 
experimental  data.  The  primary  focus  in  this  project  was  on  analyzing  data  from  the  ONR-sponsored 
SP1CE04/LOAPEX  experiment.  The  SPICE04  data  set  consists  of  measurements  from  two  vertical  line 
arrays  deployed  in  the  North  Pacific.  In  2007-2008  horizontal  line  array  data  from  SWellEx-96,  another 
ONR-sponsored  experiment,  was  also  analyzed  using  the  multitaper  method. 

In  addition  to  the  principal  investigator  (PI)  Kathleen  Wage,  three  master’s  students  at  George  Mason 
University  worked  on  this  project.  Songshun  Xu  researched  previous  work  on  single-snapshot  adaptive 
array  processing  algorithms,  and  summarized  in  a  “scholarly  paper”  (a  type  of  a  master’s  project  at 
GMU)  presented  in  summer  2008.  Khalid  AlMuhanna  explored  the  problem  of  using  ambient  noise 
measurements,  specifically  the  sample  covariance  matrix,  to  estimate  the  mode  shapes  for  a  deep  water 
waveguide.  Mr.  AlMuhanna  defended  his  master’s  thesis  in  June  2008.  Richard  Wheelock  applied 
multitaper  methods  to  analyze  the  angular  spread  of  signals  received  on  horizontal  arrays  deployed  in 
the  SWellEx-96  experiment.  Mr.  Wheelock  defended  his  MS  thesis  in  August  2008. 

WORK  COMPLETED 

SPICE04/LOAPEX  Analysis:  An  initial  analysis  of  the  250  Hz  receptions  on  the  40-element  shallow 
VLA  deployed  for  SPICE04  was  completed.  Results  of  this  analysis  were  presented  at  North  Pacific 
Acoustic  Laboratory  (NPAL)  Workshop  in  April  2006  [2]  and  at  the  fall  2006  Meeting  of  the  Acoustical 
Society  of  America  [3].  In  2008  Kathleen  Wage  and  PhD  student  Tarun  Chandrayadula  developed  a 
method  of  estimating  missing  array  navigation  data  using  an  empirical  orthogonal  function  model  of  the 
array  motion.  Mr.  Chandrayadula  presented  this  work  at  the  IEEE/MTS  Oceans  Conference  in 
September  2008  [4],  The  method  has  been  applied  to  estimate  missing  navigation  data  for  the  LOAPEX 
experiment,  facilitating  processing  of  the  acoustic  receptions. 

Multitaper  Array  Processing:  An  important  focus  of  this  project  is  an  investigation  of  the  application  of 
multitaper  spectral  estimation  methods  to  the  passive  sonar  array  processing  problem.  The  initial  work 
on  this  topic  was  summarized  in  an  invited  talk  given  at  the  Acoustical  Society  of  America  meeting  in 
June  2007  [5]  and  an  invited  paper  presented  at  the  Asilomar  Conference  on  Signals,  Systems,  and 
Computers  in  November  2007  [6].  While  the  initial  work  focused  on  applying  the  multitaper  approach 
to  equally-spaced  arrays,  this  year’s  work  investigated  how  to  apply  the  method  to  non-uniform  arrays, 
such  as  the  horizontal  array  deployed  in  the  SWellEx-96  experiment.  Master’s  student  Richard 
Wheelock  used  the  multitaper  processor  designed  for  the  SWellEx  array  to  analyze  multipath  spread  in 
the  SWellEx  data  set.  Experimental  results  were  compared  to  ray-based  simulations. 

Deep  Water  Ambient  Noise  Analysis:  The  year-long  deployment  of  long  VLA’s  as  a  part  of  the 
SPICE04  experiment  provided  an  opportunity  to  study  ambient  noise  statistics  in  the  deep  water 
channel.  This  project  analyzed  the  deep  water  noise  in  the  North  Pacific,  focusing  primarily  on  whether 
an  eigenanalysis  of  noise  covariance  matrices  revealed  anything  about  the  underlying  acoustic  modes. 
Results  indicate  that  the  eigenvectors  of  the  noise  covariance  matrix  resemble  the  low-order  acoustic 
modes.  Master’s  student  Khalid  AlMuhanna  presented  the  initial  results  of  this  work  at  the  the  IEEE 
Underwater  Acoustic  Signal  Processing  Workshop  in  October  2007  [7],  and  PI  Kathleen  Wage 
presented  the  latest  results  at  the  summer  2008  meeting  of  the  Acoustical  Society  of  America  [8]. 

Planning  for  2009  Philippine  Sea  Experiment:  In  2008  Kathleen  Wage  analyzed  the  effect  of  array 


element  spacing  on  mode  resolution  for  the  Philippine  Sea  environment.  This  work  led  to  a  proposed 
design  for  the  axial  segment  of  the  distributed  vertical  line  array  that  will  be  deployed  in  the  2009 
Engineering  Test/Pilot  study. 

RESULTS 

This  section  summarizes  results  on  the  SPICE04  analysis,  multitaper  array  processing,  and  deep  water 
noise  analysis  obtained  during  the  life  of  the  grant. 

SPICE04  Analysis:  During  the  SPICE04  experiment  a  40-element  vertical  line  array  (VLA)  received 
signals  from  two  broadband  sources  at  ranges  of  500  km  and  1000  km.  Both  sources  were  moored  at 
depths  near  the  sound  channel  axis,  thus  they  directly  excited  the  lowest  modes.  The  VLA  was  designed 
to  spatially  resolve  20  modes  at  the  center  frequency  of  250  Hz.  Figure  1(a)  shows  the  modeshapes  for 
the  environment  at  the  SPICE  VLA  location  assuming  a  frequency  of  250  Hz.  The  line  of  crosses 
indicates  the  locations  of  the  sensors  for  the  40-element  shallow  VLA  deployed  in  SPICE04.  Blue 
crosses  mark  the  good  sensors  and  red  crosses  mark  the  locations  of  the  two  bad  sensors.  To  illustrate 
the  mode-resolving  power  of  the  SPICE  VLA,  Figure  1(b)  shows  the  mode  filter  “beampattern”. 
Narrowband  mode  filters  are  commonly  designed  using  a  least  squares  criteria,  which  leads  to  the 
pseudo-inverse  (PI)  mode  filter.  For  the  PI  filter  the  vector  a (uj)  of  estimated  mode  coefficients  at 
frequency  u>  is  computed  as  follows 

a(w)  =  (ErE)-lErp  (w),  (1) 

^ - V - ' 

wr 

where  E(w)  is  the  matrix  of  modeshapes  sampled  by  the  array  and  p(u>)  is  the  measured  pressure  along 
the  array.  The  mode  filter  is  denoted  by  the  matrix  W.  In  mode  filtering,  the  beampattern  is  defined  as 
20  log10(W;  E).  The  mth  row  of  the  beampattern  matrix  corresponds  to  the  projection  of  the  modes 
into  the  estimate  for  mode  m.  The  PI  filter  is  constrained  pass  the  desired  mode  with  unity  gain  and 
place  nulls  at  the  location  of  other  modes  included  in  the  pseudo-inverse.  It  does  not  constrain  the 
output  at  modes  not  included  in  the  pseudo-inverse  [9],  The  beampattern  in  Figure  1(b)  is  for  a  20-modc 
processor,  thus  the  beampattern  shows  a  diagonal  structure  up  to  mode  20.  As  indicated  by  the  figure, 
modes  higher  than  20  are  expected  to  project  into  the  lower  order  mode  estimates.  Assuming  the 
pseudo-inverse  is  well-conditioned,  this  modal  crosstalk  is  relatively  low. 

Broadband  mode  filtering  was  implemented  for  the  SP1CE04  data  by  first  transforming  the  data, 
implementing  narrowband  PI  mode  filters  for  each  frequency  bin,  and  then  inverse  transforming. 
Corrections  for  VLA  motion  are  implemented  as  a  part  of  the  mode  filtering  process.  See  the  papers  by 
Wage  et  al.  for  more  information  about  mode  filtering  methods  [9,  10].  Figure  2  shows  an  example  of 
the  received  signals  at  the  VLA  from  each  of  the  two  sources  and  the  corresponding  estimated  time 
series  for  modes  1,  5,  and  10.  The  received  pressure  field  consists  of  a  set  of  early  arrivals  that  look  like 
up-  and  down-going  planewaves  followed  by  a  complex  finale  consisting  of  the  low  order  modes.  As 
expected,  the  reception  due  to  the  source  at  1000  km  range  shows  more  dispersive  spreading  than  the 
source  at  500  km  range.  The  finale  also  appears  to  be  more  complex  at  1000  km  range  due  to  additional 
internal  wave  scattering.  The  mode  estimates  shown  in  Figure  2  reveal  several  important  characteristics. 
First,  the  crosstalk  from  the  higher  order  modes  that  make  up  the  planewave  arrivals  is  obvious  in  each 
of  the  plots.  This  crosstalk  could  be  removed  by  time-windowing  since  it  is  well-separated  in  time  from 
the  true  low-mode  arrivals  in  the  finale.  Second,  at  the  shorter  (500  km)  range,  the  lowest  modes  are 
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Figure  1:  Modeshapes  and  the  corresponding  mode  filter  beampattern  at  the  SPICE04  VIA  location. 

The  line  of  blue  crosses  in  Figure  1(a)  indicates  the  good  sensor  positions ;  the  red  crosses  indicate 
sensor  failures.  The  mode  filter  beampattern  in  Figure  1(b)  was  computed  assuming  a  pseudo-inverse 
mode  filter  designed  for  20  modes.  For  more  information  on  mode  filter  beampatterns,  see  [9]. 


dominated  by  a  single  arrival.  At  the  longer  (1000  km)  range,  there  is  more  evidence  of  multipath 
spread  within  each  mode.  This  is  to  be  expected  since  the  longer  range  should  be  more  affected  by 
internal  wave  scattering.  Figure  3  shows  the  stacked  time  series  for  mode  1  over  a  50-day  period  of  the 
experiment.  There  is  some  obvious  time  drift  in  the  estimates  that  can  probably  be  attributed  to  motion 
of  the  sources.  (Source  motion  was  measured  and  will  be  removed  in  future  analysis.)  Figure  3 
indicates  that  mode  1  exhibits  a  relatively  sharp  cutoff  after  the  dominant  arrival  comes  in.  Note  that 
this  is  in  contrast  to  data  from  the  earlier  ATOC  experiment,  where  the  main  arrival  was  followed  by  a 
scries  of  lower  amplitude  arrivals  [10].  These  weak  late  arrivals  were  attributed  to  bathymetric 
scattering  near  the  ATOC  source  located  on  Pioneer  Seamount.  Assuming  that  bathymetric  scattering  is 
responsible,  we  would  not  expect  to  see  these  weak  late  arrivals  when  the  source  is  suspended  in  the 
water  column,  as  the  SPICE04  sources  were.  The  absence  of  these  late  arrivals  simplifies  the  analysis  of 
the  SP1CE04  data  set,  making  it  easier  to  determine  the  fluctuations  due  to  internal  waves. 
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Figure  2:  Receptions  at  the  shallow  VLA  from  the  sources  SI  and  S2  during  the  SPICE04 
experiment.  The  top  plots  in  each  column  show  the  received  pressure  as  a  function  of  time  and  depth. 

The  early  arrivals  consist  of planewaves  and  the  late-arriving  energy  consists  of  the  lowest  order 
modes.  As  expected  the  receptions  at  1000  km  range  show  more  dispersive  spreading  than  the 
receptions  at  500  km  range.  The  lower  three  plots  in  each  column  show  the  estimated  time  series  for 
modes  /,  5,  and  10.  Note  the  presence  of  crosstalk  from  the  early  planewave  arrivals  into  the  modes. 


S1-HLF5:  Mode  1  receptions 


S2-HLF5:  Mode  1  receptions 


(a)  Mode  1  arrivals  at  500  km  range. 


(b)  Mode  1  arrivals  at  1000  km  range. 


Figure  3:  Comparison  of  the  estimated  time  series  for  mode  1  for  the  SI  ( Figure  3(a))  and  S2 
(Figure  3(b))  sources  in  the  SPICE04  experiment  The  plots  show  the  receptions  over  50  days  of  the 
experiment  Note  that  source  motion  has  not  been  removed  from  these  plots,  thus  there  is  a 
noticeable  drift  of  the  arrival  times,  particularly  for  source  SI  at  500  km  range . 


Multitaper  Array  Processing:  The  multitaper  method,  formulated  by  Thomson  for  nonparametric 
spectral  estimation  [1],  is  widely  used  in  the  analysis  of  geophysical  time  series,  but  only  a  few 
researchers  have  applied  it  to  the  spatial  spectral  estimation  problem,  e.g.,  [11,  12,  13].  The  basic  idea  is 
to  apply  a  set  of  orthonormal  tapers  (windows)  to  the  data  and  then  average  the  spectra  obtained  with 
these  windows.  Averaging  over  different  tapers  reduces  the  variance  of  the  spectral  estimate,  in  the 
same  way  that  averaging  over  snapshots  in  the  conventional  Welch-Bartlett  approach  [14]  reduces  the 
variance. 

Figure  4  illustrates  a  basic  framework  for  multitaper  array  processing.  As  shown,  multitaper  processing 
can  be  viewed  as  a  series  of  operations.  First,  the  data  is  passed  through  a  beamspace  processor  that 
projects  the  data  into  several  orthogonal  beams  centered  around  the  angle  of  interest  6.  The  beamspace 
is  defined  by  a  set  of  orthogonal  tapers  such  as  the  discrete  prolate  spheroidal  sequence  (DPSS)  tapers 
suggested  by  Thomson.  As  an  example,  Figure  5  shows  the  beam  response  of  a  set  of  four  DPSS  tapers. 
To  compute  an  estimate  of  the  continuous  spectrum,  the  outputs  of  the  beamspace  processor  are 
averaged,  i.e.,  the  multitaper  estimate  is 


Smt{0)  =  Y.  Qfcl<7<:(fl)l2>  (2) 

k=i 

where  n^.  is  the  weight  given  to  the  k  beam  output.  The  weights  can  be  set  to  ^  for  a  simple  average  or 
determined  adaptively  (see  Thomson’s  seminal  paper  [1]  for  additional  details).  Equation  2  provides  an 
estimate  of  the  continuous  spectrum,  but  as  Thomson  suggests,  line  components  need  to  be  estimated 
separately.  In  the  passive  sonar  problem,  line  components  correspond  to  planewave  arrivals.  These  can 
be  estimated  using  a  linear  regression  on  the  output  of  the  beamspace  processor.  Drosopoulos  and 
Haykin  discuss  the  detection  of  line  components  in  detail  for  a  radar  problem  [11],  To  work  with  sonar 
data,  where  the  signal  is  typically  modeled  as  complex  Gaussian,  these  techniques  must  be  modified 
slightly,  as  discussed  in  [6].  An  iterative  scheme  is  useful  for  detecting  low-level  planewave  components 
in  the  presence  of  loud  interferers.  After  one  pass  through  the  beamspace  processor,  the  planewave 
signals  are  estimated  and  removed.  The  residual  signal  is  passed  through  the  beamspace  processor 
again,  and  the  planewave  detection/estimation  algorithm  is  repeated.  Figure  6  illustrates  the  output  of 
the  multitaper  processor  after  each  iteration  when  the  input  signal  consists  of  two  closely-spaced 
planewaves.  For  this  simulation  the  array  had  128  elements  with  half-wavelength  spacing.  The  two 
sources  had  different  power  levels:  the  first  had  20  dB  signal-to-noise  ratio  (SNR)  with  respect  to  a 
white  noise  floor  and  the  second  had  -10  dB  SNR  with  respect  to  the  noise  floor.  As  Figure  6  shows, 
when  no  planewave  detection  is  implemented,  the  sources  appear  as  a  blurred  peak  in  the  output 
spectrum.  After  the  first  iteration,  the  strong  source  is  detected,  but  the  weak  source  still  appears  to  be 
blurred.  After  the  second  iteration,  both  the  strong  and  weak  sources  are  clearly  visible  in  the  output. 

Figure  7  compares  the  performance  of  the  multitaper  processor  with  the  performance  of  the  standard 
minimum  power  distortionless  response  (MPDR)  processor  for  a  complex  scenario.  (See  the  text  by 
Van  Trees  for  a  thorough  description  of  the  MPDR  method  [15].)  In  this  example,  the  array  consists  of 
128  elements  with  half-wavelength  spacing.  The  environment  contains  one  moving  source  and  13 
stationary  interferers.  The  moving  source  has  a  -10  dB  SNR  with  respect  to  the  white  noise  floor;  the 
stationary  interferers  have  varying  power  levels.  Note  that  the  interferers  located  near  endfire  (four  near 
forward  endfire  and  four  near  aft  endfire)  are  perfectly  coherent.  Each  processor  uses  4  snapshots  to 
compute  its  estimates;  the  MPDR  processor  requires  diagonal  loading  to  stabilize  the  inversion  of  its 
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Figure  4:  Multitaper  processing  framework .  The  multitaper  approach  projects  the  data  into  an 
orthogonal  beamspace  and  computes  the  spectral  estimate  by  averaging  over  the  beam  outputs .  An 
iterative  algorithm  is  used  to  remove  planewave  (line)  components. 


sample  covariance  matrix.  Figure  7  shows  that  the  multitaper  processor  produces  spatial  spectra  with 
significantly  lower  variance  than  the  MPDR  processor.  As  a  result  the  low-level  moving  source  is  easier 
to  see  in  the  multitaper  output  than  it  is  in  the  MPDR  output.  In  addition,  the  multitaper  processor 
handles  the  coherent  sources  without  a  problem,  whereas  the  MPDR  processor  nulls  out  some  of  the 
coherent  signal  components.  While  forward-backward  averaging  can  help  the  MPDR  processor  with 
coherent  sources,  this  requires  additional  computation  that  the  multitaper  approach  does  not  require.  As 
Figure  7  illustrates,  the  multitaper  processor  performs  extremely  well  for  a  long  array  with  very  low 
numbers  of  snapshots.  Additional  work  is  ongoing  to  quantify  the  performance  of  the  multitaper 
processor  and  to  adapt  the  multitaper  method  to  process  non-planewave  signals. 

In  2006-2007  we  showed  that  the  multitaper  approach  works  well  with  low  numbers  of  snapshots,  and 
can  reliably  detect  low-  level  planewave  arrivals  in  a  complicated  background  containing  coherent  and 
incoherent  interference  sources  [6].  The  main  accomplishment  in  2007-2008  was  applying  the 
multitaper  method  to  analyze  real  horizontal  line  array  data  from  the  ONR-sponsored  SWellEx-96 
experiment  [16].  Since  SWellEx  used  unequally-spaced  arrays,  new  taper  designs  had  to  be 
implemented.  In  Thomson’s  original  derivation,  the  sample  rate  is  assumed  to  be  uniform,  and  the 
resulting  tapers  correspond  to  the  discrete  prolate  spheroidal  sequences  (DPSS).  Bronez  extended 
Thomson’s  spectral  estimation  method  to  irregularly  sampled  multidimensional  processes  by  deriving  a 
set  of  generalized  prolate  spheroidal  sequences  (GPSS)  [17,  18].  The  derivation  of  the  GPSS  is 
summarized  as  follows.  The  goal  in  multitaper  approach  is  to  estimate  the  integrated  spatial  spectrum, 
i.e., 

I  Mmax 

1’a  =  7T  S(kz)dkz, 

2,7V  J  Am[n 

where  the  analysis  band  A  determines  the  spatial  resolution  of  th  estimator.  Bronez  focused  on 
quadratic  estimates  of  the  power  spectrum,  i.e., 

^  K 

Pa  —  I^Pl2  =  PH  WWHp  (p  =  received  data  vector). 

k  i 

Bronez  determines  an  optimal  set  of  tapers  that  1)  guarantee  an  unbiased  estimate  when  the  spatial 
spectrum  S(k2)  is  flat  across  the  band,  2)  minimize  the  variance  of  the  spectral  estimator,  and  3) 


Spatial  response  of  tapers  1  and  2 


Figure  5:  Beam  response  of  four  discrete  prolate  spheroidal  sequence  tapers  used  in  a  typical 

multitaper  spatial  processor ; 


minimize  amount  of  “global”  bias  (due  to  signals  outside  analysis  band).  Using  these  criteria,  the 
resulting  optimization  problem  can  be  written: 

Min  ( w"Rbw*  -  Y1  w”RAwfc)  s.t.  Y1  w^lRBw*  =  j/^max  ^min-, 

\k= 1  fc=l  /  fc=l  Z7T 

Then  the  GPSS  tapers  are  a  solution  to  the  following  generalized  eigenvalue  problem: 


RaWic  =  AfcRBWfc. 


where 

RA(n,m)  —  f  e~ikl(Zn~Zm)dkz  RB(»,m)  —  (  (rik‘(Zn~tm)(lkz 

2n  J a  2ix  J b 

B  represents  the  signal  band,  i.e.,  the  band  over  which  signals  could  be  arriving.  For  equally-spaced 
arrays,  RB  ~  I,  yielding  a  standard  eigenvalue  problem  and  the  DPSS  tapers.  If  K  of  these  GPSS 
tapers  are  used  in  the  spectral  estimate,  the  resulting  bias  and  variance  of  the  estimate  are  given  below: 
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Figure  6:  Example  of  iterative  detection  of planewave  components.  The  plots  show  results  for  a 
128-element  array  with  half-wavelength  spacing.  Two planewaves  with  SNR’s  of  20  dR  and  -10  dB 
relative  to  a  white  noise  floor  are  impinging  on  the  array.  With  no  planewave  detection,  the  output 
spectrum  is  the  blurred  peak  shown  in  the  top  plot.  The  results  of  the  first  detection  iteration  (shown 
in  the  middle  plot)  indicate  that  the  strong  planewave  is  detected.  After  the  second  detection  iteration 

(bottom  plot),  both  planewaves  are  visible. 


Figures  8  and  9  illustrate  the  application  of  the  GPSS  design  method  to  the  28-element 
unequally-spaced  horizontal  line  array  (HLA)  used  in  SWellEx-96.  The  left  plot  of  Figure  8  shows  the 
generalized  eigenvalues  as  a  function  of  the  resolution  bandwidth,  and  the  right  plot  shows  the 
corresponding  sidelobe  leakage.  The  plots  show  the  results  for  the  first  6  GPSS  tapers.  As  indicated 
above,  the  bias  (or  leakage  from  outside  the  analysis  band)  goes  down  as  the  eigenvalue  A*  increases, 
thus  the  tapers  that  should  be  used  are  the  ones  whose  eigenvalues  are  close  to  1 .  For  the  S  WellEx  array 
operating  at  200  Hz,  two  good  tapers  will  be  obtained  if  the  resolution  half-bandwidth  is  chosen  to  be 
u  —  0.3  (u  is  the  directional  cosine).  Figure  9  shows  these  two  tapers  and  their  corresponding 
beampattems.  Note  that  the  two  tapers  are  essentially  non-zero  over  different  parts  of  the  array,  thus  the 
GPSS  solution  suggests  that  sub-array  processing  provides  the  best  sidelobe  performance  for  this  array. 
Figures  10  and  1 1  illustrate  the  good  performance  (as  compared  to  conventional  and  minimum  variance 
adaptive  processing)  of  the  multitaper  method  on  real  data  from  SWellEx-96. 

Deep  Water  Ambient  Noise  Analysis:  This  project  investigated  the  idea  of  estimating  the  acoustic 
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(a)  Multitaper  processor  output. 


(b)  MPDR  processor  output. 


Figure  7:  Comparison  of  the  outputs  of  the  multitaper  and  minimum  power  distortionless  response 
processors  for  a  complex  scenario  containing  one  moving  source  and  13  stationary  interferers.  The 
array  consists  of  128  elements  with  half-wavelength  spacing .  The  color  images  show  the  estimated 
spatial  spectra  for  100  different  trials ,  and  the  line  plots  at  the  bottom  of  each  figure  show  the  results 

for  the  100th  trial. 


modeshapes  using  deep  water  ambient  noise  measurements  made  during  the  SPICE04  experiment. 
Although  noise  measurements  were  not  the  primary  focus  of  SPICE04,  the  experiment  provided  a  large 
data  set  of  measurements  in  the  band  below  200  Hz  made  with  a  40-element  vertical  line  array  (VLA). 
The  basic  idea  is  that  ambient  noise  can  be  represented  as  a  sum  of  uncorrelated  acoustic  modes,  which 
is  a  reasonable  assumption  for  distant  sources,  e.g.,  see  the  work  of  Kuperman  and  Ingenito  on  noise 
modeling  [19].  In  this  case  the  eigenvectors  of  the  noise  covariance  matrix  for  a  vertical  line  array 
should  correspond  to  the  sampled  modeshapes.  Other  authors  have  investigated  using  an 
eigendecomposition  of  the  noise  covariance  to  estimate  the  mode  functions  in  shallow  water,  e.g.,  Wolf 
et  al.  [20],  Hursky  et  al.  [21],  and  Nielsen  and  Westwood  [22].  While  the  same  approach  should  work 
for  deep  water  environments,  very  few  deep  water  experiments  have  deployed  arrays  with  sufficient 
aperture  to  resolve  the  modes.  The  work  of  D’Spain  et  al.  on  estimating  acoustic  modeshapes  using 
measurements  of  earthquake  T-phase  arrivals  is  the  only  deep  water  example  of  this  method  we  have 
found  in  the  literature. 


GPSS  eigenvalues:  SWellEx  array,  freq=200  Hz 
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Figure  8:  Generalized  prolate  spheroidal  sequences  for  SWellEx-96.  The  left  plot  shows  the 
generalized  eigenvalues  as  a  function  of  the  resolution  bandwidth  for  the  28-element  SWellEx-96 
geometry  operating  at  a  frequency  of 200  Hz.  The  right  plot  shows  the  corresponding  sidelobe 
leakage  (1  —  X/J  for  each  taper.  The  dashed  line  in  the  left  plot  indicates  the  resolution  of  an 

untapered  conventional  beamformer. 


GPSS  tapers:  SWellEx  array,  freq=200  Hz 


Taper  beampatterns:  SWellEx  array,  freq=200  Hz 


Figure  9:  Generalized  prolate  spheroidal  sequence  (GPSS)  tapers  and  their  corresponding 
beampatterns.  The  left  plot  shows  the  first  two  GPSS  tapers  designed  for  the  SWellEx-96  array 
operating  at  200  Hz.  The  right  plot  compares  the  beampatterns  for  these  tapers  with  the  pattern  for 
the  conventional  (uniformly-weighted)  beamformer.  The  sidelobe  performance  of  the  GPSS  tapers  is 
significantly  better  than  for  the  uniform  weighting,  at  the  expense  of  reduced  resolution. 
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Figure  10:  SWellEx-96  experimental  setup  and  ray  theory  predictions .  The  left  plot  shows  the 
location  of  the  transmitter  and  the  South  horizontal  receiving  array.  For  this  case  the  range  is 
2.6  km.  The  right  plot  shows  the  ray  theory  predictions  for  propagation  from  the  source  at  a  depth  of 
approximately  60  m  to  the  bottom-mounted  receiving  array.  Ray  theory  predicts  three  arrivals:  two 

refracted  paths  and  one  reflected  path . 
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Figure  11:  Comparison  of  beantformer  performance  for  the  S  Well  Ex-9  6  reception  at  a  range  of 
2.6  km.  The  three  plots  show  the  results  of  conventional  beamforming ,  multitaper  beamforming ,  and 
minimum  power  distortionless  response  (MPDR)  beamforming  [15].  Each  beantformer  had  access  to 
4  snapshots  of  data;  the  MPDR  beamformer  was  diagonally-loaded  to  stabilize  the  inverse.  The  x’s 
on  each  plot  mark  the  arrival  angles  predicted  by  ray  theory.  While  the  conventional  beamformer  is 
unable  to  resolve  the  multipath  and  the  MPDR  output  is  unreliable ,  the  multitaper  processor 
successfully  localizes  two  of  the  three  predicted  arrivals. 


This  project  focused  on  the  important  signal  processing  issues  involved  in  this  type  of  mode  inversion. 
The  work  described  in  [23]  explored  the  effects  of  analysis  frequency,  array  aperture,  array  tilt,  and 
number  of  data  snapshots  using  a  set  of  simulations  designed  to  model  the  SPICE04  environment.  The 
conclusions  of  the  simulation  study  are  as  follows.  First,  as  expected,  the  aperture  of  the  array 
determines  which  acoustic  modes  can  be  extracted  from  the  eigendecomposition  of  the  noise 
covariance.  The  longer  the  aperture,  the  more  modes  can  be  successfully  extracted.  Second,  the  analysis 
frequency  also  affects  which  modes  can  be  estimated  in  this  manner.  As  frequency  increases,  the 
number  of  propagating  modes  increases.  If  the  number  of  propagating  modes  exceeds  the  number  of 
hydrophones  on  the  array,  then  the  results  are  ambiguous.  Also  the  frequency  affects  the  turning  depths 
of  the  modes.  As  frequency  decreases,  the  turning  depths  get  farther  apart,  thus  requiring  arrays  with 
greater  aperture  in  order  to  span  the  modes.  Third,  tilt  can  have  a  significant  effect  on  the  modeshape 
estimates.  While  the  array  location  can  be  tracked  using  a  transponder  network,  it  is  difficult  to  use  the 
navigation  information  to  correct  for  tilt  because  the  direction  of  the  noise  signals  is  unknown  (in 
general  noise  is  coming  from  all  directions).  Fourth,  the  number  of  snapshots  determines  the  accuracy 
of  the  noise  covariance  matrix  estimate,  thus  the  accuracy  of  the  eigenvectors  of  that  matrix.  For  the 
SPICE04  array  geometry,  the  simulation  study  indicates  that  is  is  possible  to  reliably  estimate  the 
lowest  two  modes  for  frequencies  on  the  order  of  10-20  Hz.  The  SPICE04  VLA  spanned  1400  m  of  a 
5500  m  deep  waveguide.  With  this  span,  tilts  up  to  1  degree  across  the  array  have  no  effect  on  the  mode 
estimates.  The  simulations  showed  that  100 N  snapshots  (where  N  —  number  of  hydrophones)  is 
sufficient  to  produce  a  good  estimate  of  the  eigenstructure  of  the  noise.  The  SPICE04  noise  recordings 
contain  enough  data  that  acquiring  WON  snapshots  for  the  analysis  is  not  a  problem. 

In  addition  to  simulation  studies,  this  project  applied  the  empirical  mode  analysis  technique  to  real  data 
from  the  SPICE04  experiment.  As  an  example,  Figure  12  compares  the  empirical  modes  determined  by 
eigendecomposition  with  the  true  modes  of  the  waveguide  for  yearday  455  of  the  experiment.  The  true 
modes  were  determined  by  the  standard  Prufer  mode  code  using  the  sound  speed  profile  obtained  from 
temperature  and  salinity  measurements  made  along  the  VLA.  As  expected,  the  agreement  for  modes  1 
and  2  at  the  1 1  Hz  center  frequency  is  quite  good.  Figure  12  also  shows  the  correlation  between  the 
empirical  modeshapes  and  the  true  modeshapes,  which  is  a  useful  indicator  of  how  closely  the  estimates 
agree.  Figure  13  shows  the  empirical  modeshapes  as  a  function  of  frequency  for  the  first  two  acoustic 
modes  on  yearday  455.  These  plots  demonstrate  that  good  performance  is  achieved  for  frequencies  less 
than  20  Hz.  The  results  for  many  other  yeardays  are  similar,  though  the  method  is  not  successful  for  all 
yeardays.  Based  on  the  analysis  of  the  SPICE04  data  set,  the  eigendecomposition  approach  has 
difficulties  on  days  when  the  noise  level  is  exceptionally  high  at  low  frequencies.  These  high  noise 
events  appear  to  correlate  with  high  wind  stress  (as  measured  by  satellite),  which  could  indicate  the 
presence  of  a  storm  above  the  array.  If  there  is  a  local  storm,  the  assumption  of  distant  noise  sources 
would  be  violated,  causing  the  method  to  fail.  Further  analysis  of  the  high  noise  cases  is  ongoing. 

IMPACT/APPLICATIONS 

The  Navy  has  a  vital  interest  in  being  able  to  detect  and  localize  quiet  sources  in  the  presence  of  high 
levels  of  ambient  noise  and  strong  moving  interferers,  such  as  surface  ships.  The  expected  outcomes  of 
this  research  are  new  algorithms  for  processing  non-stationary  signals  and  new  approaches  for 
incorporating  propagation  physics  into  these  algorithms.  As  the  SWellEx-96  analysis  described  in  this 
report  demonstrates,  the  multitaper  approach  is  useful  for  analyzing  multipath  arrivals  on  horizontal  line 
arrays.  In  addition,  the  SPICE04  noise  analysis  suggests  that  analyzing  noise  data  can  provide  valuable 
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Figure  12:  Comparison  of  empirical  ntodeshapes  computed  from  an  eigenvector  decomposition  of 
the  noise  sample  covariance  matrix  with  the  modes  computed  using  the  Prufer  mode  code  [ 24 j  and 
the  measured  environmental  profile.  The  left  plot  shows  the  ntodeshapes,  and  the  right  plot  shows  the 
correlation  of  the  empirical  modeshapes  with  the  actual  modeshapes.  These  plots  show  the  results  for 
a  center  frequency  of  11  Hz.  At  this  frequency  the  first  two  empirical  modes  show  excellent 
agreement  with  the  true  modes  of  the  waveguide. 
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Figure  13:  Empirical  modeshapes  for  modes  1  and  2  as  a  function  of  frequency  (solid  lines).  The 
empirical  modeshapes  were  derived  from  an  eigendecomposition  of  a  noise  covariance  matrix 
measured  on  yearday  455  during  the  SPICE04  experiment.  For  reference  the  plots  also  show  the 
actual  modes  computed  by  the  Prufer  code  [24]  for  the  measured  sound  speed  profile  (dashed  lines). 
Agreement  between  the  empirical  modes  and  the  true  modes  is  good  for  frequencies  less  than  20  Hz • 


information  about  the  deep  water  waveguides,  i.e.,  a  measurement  of  the  modeshapes.  This  project  is 
exploring  whether  similar  measurements  can  be  used  to  estimate  the  sound  speed  profile.  The  ability  to 
determine  sound  speed  or  other  characteristics  of  the  channel  without  having  to  measure  them  directly 
is  relevant  to  the  design  and  development  of  adaptive  sonar  systems.  This  project  also  has  relevance  for 
science-related  applications.  Specifically,  the  adaptive  algorithms  developed  in  this  research  is  being 
used  to  analyze  data  from  long-range  propagation  experiments  sponsored  by  ONR.  The  results  will 
provide  useful  data  for  studies  of  propagation  physics  and  tomographic  inversions. 

RELATED  PROJECTS 

This  project  is  closely  related  to  the  Navy  STTR  N04-T01 1  Phase  II  award  that  is  a  joint  project 
between  3  Phoenix,  LLC  and  Dr.  Kathleen  Wage.  The  3  Phoenix  Pi’s  are  Dr.  Russ  Jeffers  and 
Mr.  Bruce  Gallemore.  The  STTR  award  funded  the  first  part  of  Richard  Wheelock's  master’s  work.  In 
addition,  this  project  is  closely  related  to  ONR  Award  N00014-06-1-0223,  which  is  an  Ocean  Acoustics 
Graduate  Traineeship  grant  for  Tarun  Chandrayadula,  one  of  Dr.  Wage’s  Ph.D.  students. 

Mr.  Chandrayadula  is  investigating  statistical  models  for  low-frequency  acoustic  signals  propagating 
underwater  and  designing  signal  processing  techniques  to  mitigate  fluctuations  due  to  internal  waves.  In 
addition  to  the  collaboration  with  Mr.  Chandrayadula,  the  PI  also  has  collaborations  with  members  of 
the  ONR-funded  North  Pacific  Acoustic  Laboratory  (NPAL)  group.  The  NPAL  Pi's  are  Dr.  Peter 
Worcester  (Scnpps)  and  Dr.  James  Mercer  (APL-UW).  Dr.  Wage  attends  the  yearly  NPAL  workshops 
and  went  to  sea  on  the  deployment  and  recovery  cruises  for  the  2004-2005  experiment.  She  is  using 
data  collected  during  this  experiment  to  test  new  adaptive  processing  algorithms.  Dr.  Wage  is  also 
participating  in  planning  for  the  upcoming  experiment  in  the  Philippine  Sea,  scheduled  for  2009. 
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